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Abstract. Non-axisymmetric components, such as spirals and central bars, play a major role in shaping galactic discs. An 
important aspect of the disc secular evolution driven by these perturbers is the radial migration of stars. It has been suggested 
recently that migration can populate a thick-disc component from inner-disc stars with high vertical energies. Since this has 
never been demonstrated in simulations, we study in detail the effect of radial migration on the disc velocity dispersion and 
disc thickness, by separating simulated stars into migrators and non-migrators. We apply this method to three isolated barred 
Tree-SPH N-body galaxies with strong radial migration. Contraiy to expectations, we find that as stellar samples migrate, 
on the average, their velocity dispersion change (by as much as 50%) in such a way as to approximately match the non- 
migrating population at the radius at which they arrive. We show that, in fact, migrators suppress heating in parts of the disc. To 
confirm the validity of our findings, we also apply our technique to three cosmological re-simulations, which use a completely 
different simulation scheme and, remarkably, find very similar results. We believe the inability of migration to thicken discs is 
a fundamental property of internal disc evolution, irrespective of the migration mechanism at work. We explain this with the 
approximate conservation of the (average) vertical and radial actions rather than the energy. This "action mixing" can be used to 
constrain the migration rate in the Milky Way: estimates of the average vertical action in observations for different populations 
of stars should reveal flattening with radius for older groups of stars. 
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1. Introduction 

Most local galaxies are obs erved to host a thick disc component 
(e.g.. IComeron et al. 201 1 ), including our Mi lky Way (MW), 
as first proposed by iGilmore and Reidl (119831) . Compared the 
thin disc, the thick disc is composed of older, metal-poor, al- 
pha enriched stars. An important question of modern galactic 
astronomy is how thick discs came into existence. 

One possibility is that they were born thick at high redshift 
from the internal g ravitational instabilities in gas-rich, turbu- 



component in chemistry and phase-space, while in the latter 
case it would rather be seen as a gradual, continuous tran- 
sition. Observationally, whether the MW thick disc is con- 
sistent with the former or the latter is heavi ly debated (e.g, 
Schonrich and Binnevll2009l:lBovv et ai]|2012h . 



Preexisting thin discs can be he ated fast th rough 



the effect of multip l e minor mergers (|Quinn et al.l 11993 



2011), the rate 



lent, clumpy discs (B ournaud et al.L 120091) or in the turbulent 



phase assoc iated with numerous gas-rich mergers (Br ook et al 



2004 12005). They coul d also have been created through accre 



Villalobos and Helmii 120081: iDi Matteo et all 
of which decreases with decreasing redshift. Evidence for 
merger encounters can b e seen in structure in phase-space 
of MW disc stars (e.g 



2012bl: IPurcell eta 



tion of galaxy satellites (Me za et al.L 120051: 1 Abadi et all 120031) . 



4 Gvr (jpomez et al 



where thick disc stars then have an extragalactic origin. 

Another possibility is that thick discs are created through 
the heating of preexisting thin discs. This can either hap- 
pen fast, on a Gyr timescale, in an early violent epoch, or 
as a continuous process throughout the galaxy lifetime. In 
the first case, the thick disc would appear as a clear distinct 



2011 



Min chev et al. 120091: Gomez et al 



), which could last for as long as 



2012al) . The recent work by Minch ev et al 



d2012bl) presented numerical results indicating that the veloc- 
ity dispersion of the oldest MW-disc stars cannot be obtained 
without substantial merger activity at high redshift. 

The fact that galactic discs heat with time has now been 
known for over 60 years, established from both observational 
and theoretical work. To explain the observed correlation be- 



2 



I. Minchev et al.: Migration and thick discs 



twee n the ages and velocity dispersions of sola r neighborhood 



stars, ISpitzer and Schwarzschildl (1195 lL 119531) suggested that 
massive gas clouds (then undetected) scattered stars from ini- 
tially circular, into more eccentric and inclined, orbits. Giant 
mole cular clouds were thought to be th e sole s c atterin g agents 
(e.g., iMihalas and Binnevlll98ll) until iLacevi (119841) showed 
that the observed ratio of the dispersion in the direction per- 
pendicular to the Galactic plane and that toward the Galactic 
center, <x~/cr s , was too low to be consistent with the predic- 
tions from this scattering process. This resulted in the develop- 
ment of models that incorporated the heating of the stel l ar disc 



1967 



from transient spiral structur e (IBarbanis and Woltjer 
Sellwood and Carlbergj 1984|) in addition to scattering from 
molecular clouds djenkins and Binnev , ll990Hjenkinsill992h . 
Other proposed models for the heating of stars inclu de scat- 
tering by halo black hole s ( Lacev and Ostrikeri 1 19851) or dark 
clusters ( Carr and Lacev 1987b , scattering by giant molecu- 
lar cl ouds and halo black holes ( Hanninen and Flvnnl 2002 , 
2004), and v e rtical disc h eating through "popping" star clus- 
ters dKroupai 120021 120051) or through "levitation" caused by 
the 2:2 resonance between the horizontal and vertical oscilla- 
tions of stars in an inside-out for ming disc, converting their ra- 
dial actions into vertical actions (ISridhar and Touma , Il996blah . 
An investigation of the origin of the age-velocity r elation in the 
MW was presented recently bv lHouse et al.l(l201 lb . by compar- 
ing with cosmological simulations. Interestingly, a recent work 
has shown that the age-velocity relation could arise not because 
of external perturbers or the gradual heating of a thin disc, but 
because the ve locity disper s ion of newly born stars decreases 



with redshift (Fo rbes et al.l, 1201 ll) . Finally, a viable internal 



mechanism for disc heating is the interaction of mult i ple sp iral 
density waves, as proposed by Minchev and Ouillenl (120061) or 
bar and spirals ( Minchev and Famaev , 201 Oh. Multiple waves 



are now known to exist in both observation (Elmegreen et al 



1992: Rix andZaritskv. 1995: Mt 


:idtet al.. |2009i) and sim- 


ulations (e.g., iMasset and Tagger! 


19971 


Ouillenetal. 2011; 


Minchev et al. l2012al Grand et al. 


20121 


Roskar etal. 2012b, 



therefore this heating process may be ubiquitous in disc galax- 
ies. 

Interestingly, non-axisymmetric patterns also give rise to 
another important process regarding the secular evolution of 
galaxy discs: radial migration. This process is related to 
the redistribution of angular momentum of stars within the 
disc, generally resulting in an increase in disc scale-length. 
Several radial migration mechanisms have been described in 
the literature: (i) the effect of the corotation resonance of 



transient spi r al den sity waves dSellwood and Binnevl, 12002 



Roskar et al 



_ 20121) . (ii) the e ffect of the non-linear coupl ing 
between multiple spiral waves dMinchev and Quillen , 2006 ) or 
bar and spirals ((Minchev and Famaev , 20101 iMinchevetaT 



201 la : Brunetti et al. . 201 lb. and ( i ii) perturb a tions due to 



2012b . Note 



minor mergers ( Quillen et al.i 2009 ; Bird et al. 
that in case (ii) above, spirals do not need to be short-lived 
transients for efficient migration; relati vely long-lived waves 
are expected to ex i st in barred discs (IBinnev and Tremainel 
20081: lOuillen et alll201 lllMinchev et alll2012ab. but are als o 
found in non-barred simulations ( e.g. | D'Onghia et al. 2012 ). 
Recently, IComparetta and Quillen (2012) showed that, even if 



patterns are long-lived, radial migration can result from short- 
lived density peaks arising from interference among density 
waves overlapping in radius. This "disc mixing" can give rise to 
a number of ph enomena, such as flattening in radial metallicity 



gradients (e.g., Minchev et al. , 201 la ; Pilkington et all 2012 ) 



and extended stellar density profiles (e.g. . [Roskar et al. 
Sanchez-Blazquez et all 120091: iMinchev et al.[l2012ab . 



2008b 



Recently there has been a growing conviction that radial 
migration can result in a thick disc formation by bringing out 
high-velocity-dispersion stellar populations from the inner disc 
and the bulge. Su ch a scenario was used, for exa mple, in the an- 
alytical model of ISchonrich and Binnevl (12009b . where the au- 
thors managed to explain the Milky Way thick- and thin-disc 
characteristics (both chemical and kinematical) without the 
need of mergers or any discrete heating processes. Similarly, 
the inc rease of disc t hicknes s with time found in the simula- 
tion by Roskar et al.l d2008a ) has been attributed to m i gratio n 
in the works bv lSales et al.l ( 2009 ) and Loebman et al. ( 201 1 ). 
We have to point out, however, that how exactly radial migra- 
tion affects disc thickening in dynamical models has never been 
demonstrated. 

In this paper we would like to examine in detail the effect of 
radial migration on the increase of radial and vertical velocity 
dispersion and, thus, on its possible relation to the formation 
of thick discs. For this purpose we study three isolated galaxy 
simulations and three simulations in the cosmological context. 
We must note that, while the simulations we examine address 
some of the most widely used initial conditions for simulating 
galactic discs: preassembled N-body discs and discs growing 
in the cosmological context, in this paper we are not testing the 
validity of other m ode ls of formation, such a s the one used by 
Sales et al.l d2009l) and lLoebman et all d201lb . 



2. Description of the simulations 

2.1. Tree-SPH N-body simulations 

For the majority of this paper we study three main runs of iso 



lated disc galaxies from the GalMer database dDi Matteo et al 



2007b : the giant SO, Sa, and Sbc runs (hereafter gSO, gSa, and 
gSb). In GalMer, for each galaxy type, the initial halo and the 
optional bulge are modeled as Plummer spheres, with charac- 
teristic masses M# and Mb, and characteristic radii rn and rg, 
respectively. Their densities are given by 

, x-5/2 



PH,B( r ) = 



3M, 



I LB 



4nP 



1 + 



(1) 



On the other hand, the initial gaseous and stellar discs fol- 
low Miyamoto-Nagai density profiles with masses M g and M s , 
and vertical and radial scale-lengths given by b g and a g , and b s 
and a s , respectively: 



(2) 



The parameters for the initial conditions of the three runs, 
including the number of particles used for the gaseous disc, 
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Fig. 1. First column: m=2 Fourier amplitudes, Ai, as a function of time (color bar) for our model galaxies. Columns 2-3: edge- 
on view for each disc and bulge (separately) for earlier and later times of their evolution, as indicated in each panel. Columns 
4-5: scale-heights, h, at ~ 2.5 disc scale-lengths. The dotted red line shows a single exponential fit. Note that in 3-4 Gyr of 
evolution with strong perturbation from bars and spirals, h only doubles despite the efficient stellar radial migration. 



Table 1. Galaxy modeling parameters 



with larger bars giving rise to stronger mi xing (also consistent 





gso 


gSa 


gSb 


M B [2.3 x 10 9 M o ] 


10 


10 


5 


M H [2.3 x 10 9 M Q ] 


50 


50 


75 


M s [2.3 x 10 9 M o ] 


40 


40 


20 


M g /M s 




0.1 


0.2 


r B [kpc] 


2 


2 


1 


r H [kpc] 


10 


10 


12 


a s [kpc] 


4 


4 


5 


b s [kpc] 


0.5 


0.5 


0.5 


a s [kpc] 




5 


6 


b g [kpc] 




0.2 


0.2 


N s 




80000 


160000 


N, 


320000 


240000 


160000 


N DM 


160000 


160000 


160000 



N g , the stellar disc and bulge, N s , and the dark matter halo, 
Nqm, are given in Table Q] The initial Toomre parameter of 
both stars and gas is taken to be Q = 1.2 as the initial con- 
dition of the Tree-SPH simulations an d particle velociti es are 
initialized with the method described in lHernauisj(ll993l). Ful l 
det ails of the simulations are given in |Pi Matteo et al.l (120071) 
and lChilingarian et al. I (l2010l) . 

These ga l axies have been studied extensively 



in 



Minchev et al.l d20 1 1 al 120 1 2al) . In these Tree-SPH N-body sim- 
ulations the gSb model de velops a bar similar to the o ne seen 



in the Milky Way (MW) dMinchev et all I2007L l2010l) . while 



both the gSO and gSa bars are substantially larger than that (see 
Fig. Q]). We have shown ( Minchev et all |201 la ) that this has 
strong effects on the migration efficiency found in these discs, 



with studies of diffusion coefficients, e.g., Bmnetti et al. 201 1 



2012 



Shevchenkoll201 II) . Most recently, we dMinchev et al 
demonstrated that, in these same models, the strongest changes 
of angular momentum occur near the bar's corotation (CR) 
region, in contradiction to the expectat ion that patterns need 



to be transient for efficient migration dSellwood and Binnev 
20021) . 



2.2. Cosmological re-simulations 

We concentrate on the controlled models just described for 
most of the analyses we are about to do, since these discs are 
easier to analyze given their constrained nature, lack of gas ac- 
cretion and mergers. However, we also will use much more re- 
alistic simulations, following the self-consistent assembly of 
galactic discs in the cosmological context, in order to test the 
validity of o ur results. These nu merical experiments were first 
presented bv lMartig et al.l (120121) . where the authors studied the 
evolution of 33 simulated galaxies from z — 5 to z — using 
the zoom-in technique described in iMartig et al.1 (12.0091) . This 
technique consists of extracting merger and accretion histories 
(and geometry) for a given halo in a A-CDM cosmological sim- 
ulation, and then re-simulating these histories at much higher 
resolution (150 pc spatial, and 10 4 ~ 5 M Q mass resolution), re- 
placing each halo by a realistic galaxy containing gas, stars 
and dark matter. We have chosen to examine three simulations, 
which we refer to as the CI, C2, and C3 models, with approx- 
imately flat rotation curves and circular velocities V c ~ 200, 
180, and 200 km/s, respectively. The bulge/bar/disc decompo- 
sition and disc mass growth for the three galaxies are presented 
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Fig. 2. First row: Time evolution of the radial vertical dispersion profiles, cr z (r), for the gSb, gSa, and gSO models. Different 
colors and line styles correspond to the times indicated in the first row, rightmost panel. The bars' corotation (CR) and outer 
2: 1 Lindblad resonance (OLR) are shown by the dotted and solid vertical lines, respectively, with colors corresponding to the 
different times. Second row: Same but for the radial velocity dispersion, cr r (r). Third row: Same but for the ratio cr z /cr r . The 
slope in cr z /cr r outside the bar shows that the radial profile flattens more with time than the vertical one, while its overall decrease 
means that the disc heats more radially than vertically. 



in the third row of Fig. 25 (CI), the f irst row of Fig. 29 ( C2), 
and the bottom row of Fig. 27 (C3) bv lMartig etaD (120121) . 

These simulations differ in several ways from the gSb, 
gSa, and gSO models described above: (i) the discs grow self- 
consistently as the result of cosmological gas accretion from 
filaments and gas-rich mergers, as well as merger debris, (ii) 
the simulation techniques are different both in the way grav- 
itational forces are computed (Particle Mesh vs hierarchical 
tree) and in the way gas dynamics is treated (sticky-particles 
vs SPH), and (iii) the spiral and bar instabilities are present for 
a cosmological time, thus the effects we see are not the result 
of recent bar formation, for example. 

3. Disc thickening and velocity dispersions 

Up until Section [8] we only analyze our constrained, Tree- 
SPH N-b ody simulations describ ed in Section l2~T1 As we have 
shown in Minchev et al. I d201 2ah . all of these models exh bit 
multiple patterns, giving rise to strong radial migration. To 
show the strength of the discs' asymmetric components, in the 
first column of Fig. [T]we plot the m = 2 Fourier amplitudes, 
A2/A0, as a function of radius for each galaxy, where Aq is the 
axisymmetric component and m is the multiplicity of the pat- 
tern. Different color curves show the time evolution during the 
first Gyr. Ao indicate the bi-symmetric structure in the disc. For 
example, for the gSO model the bar is identifiable by the smooth 



curve in the inner disc, which at later times has a maximum at 
~ 3 kpc and drops almost to zero at ~ 8 kpc. Deviations from 
zero seen beyond that radius are due to the spiral structure. 
Note that the spirals are strongest for the gSa model, related 
to the initial 10% gas fraction in its disc. The gSb model is the 
one with parameters closest to the MW: bar size ~ 3 - 4 kpc, 
v c ~ 220 km/s, and a small bulge. The circular velocities and 
ba r pattern speeds for all three models can be found in Fig. 2 
bv lMinchev et al. I d2012ah . 

To see how much discs thicken, in columns 2-3 of Fig.Q]we 
plot the edge-on view for each disc and bulge (separately) for 
earlier (0.25 Gyr) and later (4 Gyr for gSb and gSO, and 3 Gyr 
for gSa) times of their evolution, as indicated in each panel. 
Significant thickening is seen for all models, especially for gSa 
and gSO due to their larger bars. The bulges are also found to 
expand, not significantly so for the gSb. Although from these 
plots it is evident that the discs thicken, it is important to see 
how this is reflected in their scale-heights. These are shown 
in columns 4-5 at a radius 12 kpc (~ 2.5 disc scale-lengths), 
which, when scaled, is comparable to the solar distance from 
the Galactic center. The dotted red line shows single exponen- 
tial fits. 

Surprisingly, we find that, despite the strong radial migra- 
tion, the initial scale-heights, h « 0.3 kpc, hardly double in the 
3-4 Gyr of evolution for gSa and gSO, respectively, and increase 
by only ~ 50% in 4 Gyr for gSb. 
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gSb,0.5<t<1 .5 Gyr 0.5<t<4.0 Gyr 

r=12, Ar=2.00 kpc 




5 10 15 

r InIKol t k P c ] 



gSo,0.5<t<1 .5 Gyr 
r\=12, Ar = 2.00_kpi 
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0.5<t<3.0 Gyr 




5 10 15 
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Fig. 3. First row: Initial (dashed black), final (solid blue), and net (dash-dotted red) vertical velocity dispersions, cr z j, cr z j, and 
Act., respectively, as functions of the initial radius, ro, for particles ending up in the green bin at the final time indicated in each 
panel. Two time periods for each of the three models are shown, as indicated. Migrating stellar samples decrease or increase their 
velocity dispersion is such a way as to match the value at the destination. Second row: Same as the first row, but for the radial 
velocity dispersion. Third row: The corresponding initial radial distributions of stars ending up in the green bin at the final time 
indicated in each panel. The bar's CR and OLR, measured at the final times, are shown by the vertical dotted and solid black 
lines. Most migrators come from the inner disc, with a strong contribution from particles originating near the bar's CR. 



The disc thickening we see in Fig. Q] must be related to 
an increase in the vertical velocity dispersion. To see how our 
discs heat, in the first two rows of Fig. [2] we show the time 
evolution of the vertical, <r z (r), and radial, cr,.(r), velocity dis- 
persion profiles for the stellar populations of the gSb, gSa, and 
gSO models in the first, second and third columns, respectively. 
Initial bulge components are not considered. In each panel 
the different color curves show different times, as indicated in 
the top left panel. The dotted and solid vertical lines indicate 
the radial positions of the bar's corotation (CR) and 2:1 outer 
Lindblad resonance (OLR) with colors matching the different 
times. 

An interesting observation is that the <x,. profiles flatten at 
later times, especially for the gSO and gSa models, i.e., the 
galaxies with the larger bars and stronger migration histories. 
This is reflected in the ratio cr z ja r , shown in the bottom row, 
where the later curves (blue and black) develop a negative 
slope. Moreover, cr r increases more in the outer parts of the 
disc, compared to the region near the bar's CR. In fact, we 
observe radial cooling near the CR for both the gSO and gSa 



models: the solid-black curve falls below the dashed-green and 
dash-dotted blue ones. The reason for this is revealed in the 
next Sections. 

4. Cooling/heating of stellar samples migrating 
outwards/inwards 

Stellar samples in the inner galactic discs possess high velocity 
dispersions in all three components. It has been expected that 
by migrating radially outwards, stars retain their random ener- 
gies and thus thicken the disc. If this were true, then how do 
we explain the lack of sufficient disc thickening even for the 
strong large-scale disc mixing seen in Fig. [Tf To answer this 
question, we consider the samples of stars used to estimate the 
scale-heights in Fig. Q] In the first row of Fig. [3] we plot the 
initial (dashed black), final (solid blue), and net (dot-dashed 
red) vertical velocity dispersions, o% ,, cr z j, and Act-, respec- 
tively, as functions of the initial radius, ro, for particles ending 
up in the green bin (1 1 < r < 13 kpc) at the end of the given 
time interval. This is done for each galaxy in the time intervals 
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0.5 < t < 1.5 and 0.5 < t < 3 Gyr (0.5 < t < 4 Gyr for 
gSO and gSb), as indicated in the figure. The second row shows 
similar plots but for the radial velocity dispersion. The distri- 
bution of initial radii of stars ending up in the green bins at the 
corresponding times are shown in the third row. 

From the first two rows it is immediately apparent that stel- 
lar samples arriving from the inner discs decrease their velocity 
dispersions by as much as 50% when coming from the small- 
est radii, with the effect diminishing the smaller the radial dis- 
tance traveled. Conversely, stars coming from radii exterior to 
the green bin increase their velocity dispersion. Therefore, on 
the average, newcomers change their random energies in such 
a way as to approximately match the non-migrating population 
at the radius at which they arrive. Thus, even though at later 
times (second, fourth and sixth columns) a larger fraction of 
migrators enters the annulus under consideration (as seen in 
the bottom row), this has no effect on the velocity dispersion 
increase. 

The effect we just described is more dramatic for the ra- 
dial velocity dispersion: for the gSO (dissipationless) model the 
non-migrating stars (blue line portion inside the green bin) are 
hotter at the final time than the migrators, as evident by the pos- 
itive slope of the solid-blue curve. Since the majority of stars 
arrive from the inner disc, the final annulus is being cooled, 
rather than heated by the migrators! We can see that, for all 
three models, stars arriving from inside the bar's CR are radi- 
ally cooler than the non-migrators: note the positive slope in 
the blue curve inside the CR region (i.e., to the left of the dot- 
ted vertical line). Because for gSO the radial bin considered 
happens to be between the CR and the OLR, we see the most 
cooling for this model. We find similar results when the bin is 
placed at a similar location with respect to the bar's resonances 
in the case of the gSa and gSb models (see Fig. [6). 

The velocity distribution of migrating stars before and af- 
ter they have migrated are smooth and appear approximately 
gaussian, as can be seen in Fig.|9]in the Appendix. 

5. Separating migrators from non-migrators 

To assess the contribution of migrating stars to the intrinsic ve- 
locity dispersion of the disc we devise a simple procedure for 
separating migrators and non-migrators. 

(1) We consider a disc annulus of a certain width at some initial 
and then final time. 

(2) We separate all stars into "non-migrators" (those that are 
found in the selected radial bin at both the initial and final 
times) and "migrators" (those that were not present in the bin 
initially but are there at the final time). 

(3) We do this for radial bins over the entire extent of the disc, 
overlapping every 0.3 kpc, for example, in order to get infor- 
mation at every radius. 

In step (1) above we need to adopt a certain width for each 
annulus considered. We chose this as follows. Epicycle sizes 
are on the order of a ~ cr r /K, where k is the epicyclic frequency. 
In order to contain non-migrating stars, we allow for particles 
to oscillate around their guiding radii by +a, therefore, we se- 
lect our bin size to be 2a. Since both cr r and k are functions of 
radius, so is the value of a. In Fig.[10]we plot all three functions 



versus galactic radius. Note that the bin sizes increase with ra- 
dius, where the effect is stronger for the hotter gSa and gSO. 

We apply the above procedure to our gSb and gSa mod- 
els and display the results in Fig. [4] The first row shows con- 
tour density plots of the final versus initial radii for all parti- 
cles in different time spans for the gSb or gSa models, as in- 
dicated. We exclude a small portion of the inner disc (white 
rectangle in the lower left corner) in order to display better the 
outer contours. The orange dashed line shows the locus of non- 
migrating, cold particles. Nonuniform deviations from this line 
indicate that there exist preferential radii where migration ef- 
ficiency is stronger, however, these are not clearly seen. We 
note that the particle distributions presented here are not the 
true density, but a combination from all migrators and non- 
migrators resulting from the overlap of the radial bins we con- 
sider. However, deviation from the true density are small, as 
shown in Fig.QT]and discussed in the Appendix. 

The second row of Fig. [4] presents the non-migrating stars 
extracted from the distributions shown in the first row as out- 
lined in (1), (2), and (3) above. Deviation from the orange line 
are about 1-2 kpc, increasing outwards, in agreement with the 
bins-size variation shown in Fig.fTOl 

The third row of Fig.|4]shows the migrators only, obtained 
by subtracting the distributions shown in the second row from 
those in the first row. The dotted-red and solid-blue vertical and 
horizontal lines indicate the location of the bar's CR and OLR 
at the initial (horizontal) and final (vertical) times. We can now 
see clearly overdensities on each side of the orange dashed line. 
We show how our results change by considering different bin 
size values in the Appendix (see Fig . [T3jl . 

In the fourth row of Fig. |4]we plot number density contours 
of the changes in the specific angular momentum, AL, ver- 
sus the initial specific angular momentum, L, estimated during 
each time periods. Both axes are divided by the rotational ve- 
locity at each radius, therefore L is approximately equal to the 
initial radius and AL gives the distance by which guiding radii 
change. The broken-red and solid-blue vertical strips indicate 
the location of the bar's CR and OLR, where their widths cor- 
responding to the change during each time period related to the 
bar's slowing down. For all time periods considered we find a 
nice correlation with structure across the bar's CR, where stars 
inside the CR migrate outward and those outside it migrate in- 
ward. The pink lines connect some corresponding groups of 
migrated stars in the r/„ M / - r ,,„■„•„/ and the L - AL planes. At 
the earlier times shown (0.5 < t < 1.5 Gyr, left column) the 
changes in angular momentum are the strongest due to the po- 
tent spirals and their interaction with the bar. At the intermedi- 
ate times, as the spirals weaken, the typical feature of negative 
slope across the bar's CR is much better seen, in addition to 
some outer structure related to the spirals. Finally, at the later 
times (3 < t < 4 Gyr, third column) the bar's CR is clearly 
dominating. In general, although here we only show the en- 
tire simulation time for the gSb model and only the beginning 
of gSa, clumps across the bars' CR are invariably seen for all 
time periods and a ll three models. This co ntinuous bar activity 
was shown also bv lMinchev et al.l (I2012al) (see their Fig. 7). 

Contours are denser above the orange line in the r/,-, M / - 
rinitiai plane. This is most extreme during the early time period 
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Fig. 4. First row: Contour density plots of the final versus initial radii for all particles in different time spans for the gSb 
or gSa models, as indicated. The orange dashed line shows the locus of non-migrating, cold particles. Nonuniform deviations 
from this line indicate that there exist preferential radii where migration efficiency is stronger, however, these are not clearly seen. 
Second row: Non-migrating stars, extracted from the distributions above as described in the text. Third row: The migrators only, 
obtained by subtracting the distributions in the second from those in the first row. The dotted-red and solid-blue vertical/horizontal 
lines in all panels indicate the positions of the bar's CR and OLR, respectively, at the initial/final times. Strong clumps can now 
be seen across the bar's CR for example. Fourth row: Number density plots of the changes in angular momentum, AL, versus 
the initial angular momentum (or equivalently initial radius), L, divided by the rotational velocity at each radius. The broken red 
and blue vertical strips indicate the location of the bar's CR and OLR, where their widths corresponding to the change during 
each time period. The pink lines connect some corresponding groups of migrated stars. The bar's CR is active throughout the 
entire simulations. 



shown for the gSa model (rightmost column of Fig. @). The 
reason for this behavior is that, in the outer disc, stars are pref- 
erentially shifted outwards due to the decreasing stellar density. 
The effect is much stronger for gSa because the outward trans- 
fer of angular momentum fo r this simulation is muc h stronger 
than for the gSb disc, which iMinchev et al.l d2012al) related to 
the very efficient non-linear coupling among the bar and spiral 
waves of different multiplicity. Note that, since the total an- 



gular momentum of the disc must be conserved (with a small 
contribution to the halo and bulge), the outward transfer of an- 
gular moment needs to be balanced. This balance is achieved 
by transferring large amounts of mass inward of the bar's CR 
(with small negative changes in stellar guiding radii), in con- 
trast to the exponentially decreasing density in the outer disc 
(but much larger increase in guiding radii). An indication of 
this process may be found in the strong clump of inward mi- 
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Fig. 5. Migrators' contribution to the disc heating. First row: 
radial profiles of the vertical velocity dispersion, cr z {r), for non- 
migrators (dashed black curve), outward going migrators (solid 
brown curve), and inward migrators (dash-dotted green curve) 
for the gSb model, estimated in the time periods indicated in 
the top row of each column. Second row: cr z (r) estimated for 
non-migrators (dashed black curve), all migrators (i.e., inward 
and outward going, solid green curve), and the total population 
(dashed-dotted blue curve). Third and Fourth rows: same as 
the first and second rows, respectively, but for the radial veloc- 
ity dispersion profile o>(r). An unexpected result is that out- 
ward migrators cool the disc inside the bar's OLR (and inward 
ones heat it), unlike in the case of cr z . Some increase in cr : and 
<r r from migrators is seen beyond the OLR, while the inner disc 
is cooled both radially and vertically. Fig.|6]quantifies better the 
migrators' contribution to cr_ and cr r . 



grators just outside the bars' CR, seen in all plots but the latest 
time period for gSb; the reason for the latter is due to the sat- 
uration of migration efficiency in the outer disc in the absence 
of gas accretion (and thus spiral strength). 

6. Migrators' contribution to the disc velocity 
dispersion 

By separating stars in the entire disc into "migrators" and "non- 
migrators", we are now in a position to investigate their indi- 
vidual effects on the disc velocity dispersion increase. We con- 
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Fig. 6. Fractional changes in o~ z (r) and cr r (r) resulting from 
the migrated stars. The dashed black curve plots the quantity 



Acr z = (cr z all - cr. 



z,non_mig 



)/o- z>aU , where cr z>all and cr z , 



non_mig 



are 



the vertical velocity dispersions obtained from the total popula- 
tion and the non-migrators, respectively. A similar expression 
estimates Acr r , shown by the solid orange curve. The dotted-red 
and solid-blue vertical lines denote the bar's CR and OLR, cor- 
responding to the final time shown in each time period. Note 
the inverse correlation between Acr z and Acr r , indicating cou- 
pling between the radial and vertical motions of stars. 



sider the time intervals 0.5 < t < 1.5 and 1.5 < t < 3 Gyr, 
corresponding to the first two columns in Fig. |4] Because the 
heating in both radial and vertical directions mostly saturates 
around 1-2 Gyr (see Fig. the effect at later times is similar 
to what we find during the period 1 .5 < t < 3 Gyr. 

In the first row of Fig. [5]we show the radial profile of the 
vertical velocity dispersion, cr z (r) for non-migrators (dashed 
black curve), outward going migrators (solid brown curve), and 
inward migrators (dash-dotted green curve) for the gSb model. 
The dotted-red and solid-blue vertical lines denote the bar's CR 
and OLR, respectively, at the final time of the time period con- 
sidered. 

We find that inward migrators cool the disc and the outward 
ones heat it. This is consistent with the information contained 
in Fig. [3] where we found that particles migrating to a final an- 
nulus at r = 12 kpc arrive with cr : values slightly lower/higher 
than the non-migrating population (evident by the blue line's 
slightly negative slope in Fig. [3]). 



I. Minchev et al.: Migration and thick discs 



9 




5 10 15 5 10 15 5 10 15 

r [kpc] r [kpc] r [kpc] 



Fig. 7. Scale-height radial profiles for migrators, non-migrators and all stars in the time intervals 0.5 < t < 1.5 Gyr (top) and 
1.5 < t < 3 Gyr (bottom). The green solid curves indicate the fraction of migrators as functions of radius during each time period. 
The scale-height for the non-migrating stars is approximately constant with radius for all models (blue dashed curves, top). Due 
to the inner disc cooling and outer disc heating in the vertical direction from migrators (dotted-red curve, see also Fig. [6j, all 
discs flare (solid black curve). 



From the plots we just described (first row of Fig. [5} it 
is not clear how inward and outward migrators contribute to 
the velocity dispersion profile of the total population since we 
did not normalize to the number of particles in each group. 
How much cr,(r) of the total population is affected depends 
on the fraction of each kind of migrators at each radial bin. 
To demonstrate the true contribution of migrators, in the sec- 
ond row of Fig. |5] we plot the vertical velocity dispersion for 
non-migrators (dashed black curve), all migrators (i.e., inward 
and outward going, solid green curve), and the total population 
(dashed-dotted blue curve). Interestingly, we find that the ve- 
locity dispersion from both types of migrators is very similar 
to that of the non-migrating particles, except for the inner parts 
of the disc for the earlier times shown. It is surprising that not 
only do migrators not heat most of the disc, but they even cool 
it inside the CR. Some increase in o~ z can be seen outside the 
OLR. 

We now examine the effect of migrators on the radial ve- 
locity dispersion profile, <x,.(r). The third and fourth rows of 
Fig. [5] are similar to the first two, except the plotted quantity 
is o~ r (r). A distinct difference compared to the vertical velocity 
dispersion is the oscillatory behavior of <x,.(r) found in all three 
curves in the third row. Comparing the third and fourth rows for 
the earlier time period, it becomes clear that the dominant effect 
inside the bar's OLR comes from the outward migrators, which 
here cool, rather than heat the disc as in the case of cr-(r). Just 
inside and outside the OLR, the effect of migrators is to reduce 
the wiggles seen in the non-migrators' o>, consequently result- 
ing in a smooth radial velocity dispersion profile. 

To better quantify the contribution of migrators, we now 
calculate the fractional changes in cr z (r) and o>(r) resulting 
from the migrated stars for the same time periods shown in 



Figs. For the vertical velocity dispersion, we estimate these 
as 

Act z = (cr z<aU - cr z , n on_mig)/ & z,all» (3) 

where cr za ii and cr z non _ m i g are the vertical velocity dispersions 
for the total population and the non-migrators, respectively. An 
equivalent expression estimates the contribution to the radial 
velocity dispersion. We plot the results in Fig. [6] The dotted- 
red and solid-blue vertical lines denote the bar's CR and OLR, 
corresponding to the final time considered in each panel. It now 
becomes much clearer how exactly migrators affect the discs. 

For all models shown A<x z (dashed black curve) shows os- 
cillatory behavior, changing abruptly from negative inside the 
CR, to positive outside it. The maximum deviations from zero 
for gSb, for example, are ~ 7% cooling in the inner disc and 
~ 4% heating at larger radii, and appear consistent with Fig. [5] 
Clearly, such contribution to the vertical velocity dispersion 
would results in disc flaring, as we show in the next Section. 

Now we look at the migrators' contribution to the radial 
velocity dispersion, Acr r (solid orange curve). Remarkably, 
there appears to be an inverse correlation between A<x z (r) and 
Acr r (r), where the two functions cross near the bars' CR and 
OLR. This oscillatory behavior with zeros near the CR and the 
OLR is present for all three models during the first time period 
shown, which is when the discs heat the most. Note that the 
crossing point near the OLR occurs above zero as the OLR is 
shifted more outward in the disc (top to bottom), due to the de- 
creasing number of non-migrators and prevailing influence of 
outward migrators. Outside the OLR both A<x z and A<x r are pos- 
itive, with a contribution of less then 5 % (gSb) or less than 10% 
(gSa). At the later times the general trend of an inverse correla- 
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tion for the radial and vertical velocity dispersions is still seen, 
although the contribution to the Act, is larger. 

The inverse functional behavior just described indicates 
that the vertical and horizontal motions of migrators are cou- 
pled, as expected for a system governed by non-linear dynam- 
ics. 



7. The effect on the discs' scale-height 

We showed in Fig. [6] that the migrators' contribution to the 
vertical velocity dispersion increase is small and mostly in the 
outer disc. This suggests that our discs flare. To find out if this 
is really the case, we now measure the scale-height at different 
radial bins in the gSO, gSa and gSb discs. 

In Fig.|7]we show radial profiles of the discs' scale-height 
for migrators, non-migrators and all stars in the time intervals 
0.5 < t < 1.5 Gyr (top) and 1.5 < t < 3 Gyr (bottom). As in the 
figures so far, we take the initial time at 0.5 Gyr (or after the 
bars have formed) in order to avoid the effect of the initial for- 
mation of asymmetric structure. The green solid curves indicate 
the fraction of migrators as functions of radius during each time 
period. Note that during the later time period (1 .5 < t < 3 Gyr) 
the stars we refer to as "non-migrators" may have experienced 
some migration previously. This is also true for the earlier 
times, i.e., some migration takes place before 0.5 Gyr, but the 
strongest effect occurs in the range 0.5 < t < 1.5 Gyr. All our 
procedure ensures is that they stayed at the same radii during 
the later time interval. 

Fig. |7] shows that the scale-height for the non-migrating 
stars is approximately constant with radius for all models (blue 
dashed curves, top) during the earlier time period. In fact, 
for gSO the scale-height for non-migrators shows a decrease 
with radius. However, due to the inner disc cooling and outer 
disc heating in the vertical direction from migrators (dotted- 
red curve, see also Fig. |6), all discs flare (solid black curve). 
This is entirely consistent with Fig. [6] where we found that the 
contribution of migrators to the vertical velocity dispersion is 
negative inside the bar's CR and positive in the outer disc. 

By comparing the variation of the migration fraction with 
radius (green solid curve) to the scale-height of the total pop- 
ulation (black solid curve) we note an inverse correlation in- 
side the bars' CR: the more migrators present at a given ra- 
dius, the smaller the scale-height. This is especially obvious at 
the earlier time period and for the more efficient gSa and gSO. 
At larger radii, where all stars are in fact migrators (although 
fraction is smaller that 100%, these stars have migrated out be- 
fore t=0.5 Gyr), scale-height increases with migration fraction, 
which is expected. 



8. Application to cosmological re-simulations 

We now apply our newly developed technique for separating 
migrators and non-migrators to our three, much more realistic 
simulations in the cosmological context, described in Section[2] 
To assess the effect of the internal perturbers and be able to 
compare to our previous results, here we only consider the last 
Gyr of evolution before z — 0, thus, avoiding strong merger 



activity. We defer a more detailed study of these models to a 
future work. 

During the ~ 1 Gyr of evolution we study here, disc pene- 
tration by satellites occurs only at (r, f) ss (18 kpc, 13.15 Gyr) 
for CI and (r,t) * (10 kpc, 13.4 Gyr) for C3. There are no 
disc encounters within r < 30 kpc for the C2 model. All of the 
aforementioned companions have mass ratios less than 1:100. 
Face-on and edge-on stellar density contours are plotted in the 
first and second rows of Fig. [8] respectively, depicting the time 
t = 13.4 Gyr. This is roughly in the middle of the timespan 
considered for the rest of the plots in the figure. Some small 
satellites are still seen in all three simulations. For CI and C2 
these are orbiting outside 20 kpc radius, but for C3 there is a 
small companion at r ~ 10 kpc, located close to the plane at 
this time. 

8.1. First results for migration in the cosmological 
context 

To see how much migration has taken place in the last Gyr of 
evolution, in the third row of Fig. [8] for each model, we plot the 
changes in the specific angular momentum of stellar particles 
in the range 4 < r < 23 kpc and \z\ < 4 kpc. The horizon- 
tal axis plots the initial angular momentum, L, (estimated at 
t = 12.73 Gyr); this is divided by the circular velocity and, 
thus, approximately equals galactic radius. The vertical axis 
plots AL = Lf — Li, where, Lf is the angular momentum es- 
timated as t = 13.7 Gyr. Stars at AL ~ 4, for example, gain 
angular momentum which brings them ~ 4 kpc outward in the 
disc and conversely for negative AL. The dotted-red and solid- 
blue vertical lines show the radial location of each bar's CR and 
OLR. The high peak near the OLR for the C3 model is caused 
by the satellite closest to the galactic center seen in the disc 
face-on view above. 

We note that this is the first time such plots are presented 
for simulations in the cosmological context. It is interesting that 
the strongest changes in angular momentum happen nea r the 
bar's CR, just as noted recently bv lMinchev et al.l d2012al) (us- 
ing the gSb, gSa, and gSO models), despite the fact that here 
the bars have been formed for a long time (more than 6 Gyr) 
and are not transient. The strength of the bar correlates with the 
effect on the r — AL plane; the bar length is related to the pattern 
speed, as expected, with longer bars rotating slower, since bars 
exist within their CR radii. 

8.2. Contribution of migrators to the increase in 
velocity dispersion 

Now that we have established that migration is efficient in these 
runs (in fact, very si milar to our controlled simulations, see 
Minchev et al.ll2012.al Fig. 7 and Fig. [4] in this paper), we can 
look for the effect of the migrating and non-migrating stars 
(during the time period considered) on the vertical and radial 
velocity dispersion profiles, just as we did in Figs. [5] and [6] In 
the fourth row of Fig.|8]we plot the vertical velocity dispersion 
of inward (moving in) and outward (moving out) migrators, and 
the non-migrators, similarly to Fig. [5] Cooling from inward mi- 
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Fig. 8. Application to three high-resolution cosmological re-simulations. First and second rows: Face-on and edge-on density 
contours for the CI, C2, and C3 models. Third row: Change in angular momentum as a function of disc radius: effective migra- 
tion is seen for all models. Fourth row: Vertical velocity dispersion profiles for non-migrators, outward and inward migrators. 
Fifth row: Same but for the radial velocity dispersion. Bottom row: Fractional contribution of migrators to the total vertical and 
radial velocity dispersions (as Fig. [6). Results are very similar to the simpler models presented earlier. 



grators and heating from outward ones is seen, mostly inside 
the bar's CR (dotted-red vertical). This is in agreement with 
what we found for our preassembled models, with the differ- 
ence that the effect here is weaker in the outer disc, except for 
the C3 model, where a significant peak is seen between the CR 
and the OLR; the latter is most likely due to the small satellite 
mentioned earlier, taking place in the same radial range. 

The fifth row of Fig. [8]is the equivalent of the top one, but 
for the radial velocity dispersion. As in our controlled models, 



cooling from outward migrators and heating from inward ones 
(the reverse of what happens for o%) occurs between the bar's 
CR and its OLR. Even in the case of C3, which is perturbed by 
the satellite at that radius, we see the same behavior. 

Finally, the bottom row of Fig. [8] plots the fractional con- 
tribution to the vertical and radial velocity dispersions from all 
migrators estimated from Eq. [3] Similarly to Fig. |6] migrators 
contribute to A<x.(r) by ~ 5%, cooling the disc inside the CR 
for CI, just as we found for the controlled models. For C2 and 
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C3, the radial variation of Act- is more erratic, yet contribution 
is mostly less than 3%. The CI model's dynamics is dominated 
by the presence of its large bar, which could explain its closer 
similarity to the controlled simulations. 

While outward migrators for all three models are radially 
cooler than the non-migrators in the region between the CR and 
the OLR (fifth row of Fig. [8]l just as the case of our constrained 
simulations, A<r r (r) is only slightly below zero for CI and C2 
and slightly positive for C3. These differences are most likely 
due to tidal effects associated with the orbiting satellites seen 
in the face-on plots in the first row. An inverse relation between 
A<x-.(r) and A<x,.(r) for our cosmological simulations can still be 
seen. 

It is remarkable that both our isolated discs and the just de- 
scribed much more sophisticated simulations (which also use 
a completely different simulation technique) exhibit very simi- 
lar behavior. This makes our result, that migrators do not con- 
tribute significantly to the increase of the velocity dispersion of 
the disc, much stronger. We, therefore, suggest that our finding 
is a fundamental property of the secular evolution of galactic 
discs. It is possible that when stronger satellite perturbations 
are included, as is the case at the early stages of galaxy forma- 
tion, the situation may be different. 



9. Conservation of vertical and radial actions? 

Our findings so far suggest the existence of a conserved dynam- 
ical quantity for stars moving around the galactic disc, resulting 
in the conspiracy to approximately match (on the average) the 
kinematics expected at their destination. Since from Figs. 13151 
and [6] it is clear that this is not the energy, the obvious candi- 
dates are the radial and vertical actions. 



9.1. The vertical action 



Binnev and McMillan! 1 201 1 ) have shown that, to approximate 
how the vertical motion of disc stars is affected by their ra- 
dial motion in an axisymmetric potential, the adiabatic ap- 
proximation is valid, i.e. that the vertical action is adiabati- 
cally conserved in the epicyclic sens^U- As first proposed by 



Minchev et al. ( 201 lbl) and recently confirmed by Solwav et al 
d2012h . the vertical actions of migrating stars are closely con- 
served, as well. 

Note that while in an axisymmetric disc the energy and an- 
gular momentum are conserved quantities (thus, conservation 
of actions is expected), since stars migrate due to resonant per- 
turbation associated with disc asymmetries (such as the bar and 
spirals waves), the above mentioned classical integrals of mo- 
tion iarejionie£e«sarily__applicable here. Nevertheless, follow- 
ing |MjnchexetiLU (1201 lbl) . we now show that even in the case 
of a simple galactic disc model, assuming conservation of ver- 
tical action is in good agreement with our results. 



1 This means that the vertical action calculated from the vertical 
potential v P(z, t) = <t>(R(t), z) - <t>(R(t), z = 0) is very nearly conserved 
if R(t) is the horizontal motion calculated from the planar effective 
potential $>(R,0) + L:/R 2 . 



We assume that for stars at a given radius, the vertical ac- 
tion is given by J z = EJv, where E z and v are the vertical 
energy and frequency, respectively. We also assume that this 
quantity is conceived (on average) as stellar samples shift guid- 
ing radii (i.e., migrate). This is a good approximation in the 
limit that (i) the vertical motion decouples from the planar 
motion (a good assumption when z max <k r, which holds for 
most disc stars, but note that we did find evidence to the op- 
posite in Sec. [5J, and (ii) the stars move outwards much more 
slowly than their vertical oscillations (not really the case). From 
Gauss' law and Poisson's equations, v ~ ^JlnGL, therefore, 
as stars migrate radially, their vertical frequency changes as 
Vmig(f) ~ ex P( _r /2^), given that the stellar density varies as 
S ~ exp(-r/r</). To conserve the mean vertical action of a set 
of stars migrating outwards, we then require that their mean 
vertical energy, (-E z , m i g ) ~ cr\ mig , decreases as v. Therefore, the 
vertical velocity dispersion of a stellar sample migrating radi- 
ally outwards will decrease as cr z m i g (r) ~ exp(-r/4r^), where 
r is the radius at witch it ends up. 

On the other hand, the vertical velocity dispersion in the 
ambient (non-migrating) disc com ponent changes with radius 
as o~z (r) ~ VS(r)/i (e.g., eq. 7 by Ivan der Kruit and Freeman 
201 lb . Assuming the scale-height, h, remains constant, we find 
that cr z (r) ~ exp(-r/2r^). Thus, while a star's vertical energy 
does decay as it migrates outwards, it does not decay as fast as 
the underlying typical vertical energy of disc stars and, there- 
fore, still increases the vertical temperature of the disc where it 
arrives. Similarly, stars migrating inwards would cool the disc. 
This is exactly what we see in Fig. [5] at radii outside the bars' 
CR: the disc is heated/cooled by outward/inward migrators. We 
found that inside the CR the discs are, in fact, cooled verti- 
cally, since contribution is made only by inward migrators (see 
Figs.|5]and|6]l. Similarly, in the disc outskirts outward migrators 
dominate, thus the vertical temperature is increased. 

9.2. The radial action 

For the simple galaxy disc model considered above, the same 
argument which lead us to conclude that the vertical action is 
conserved can also be applied in the plane. We assume a con- 
servation of radial epicyclic action J p - E p /k where E p is the 
radial epicyclic energy and and k ~ 1 jr is the radial epicyclic 
frequency, valid for a flat rotation curve. We, therefore, expect 
that orbits conserving J p lose energy E p at the rate of 1 jr. 

There are some differences when compared to the vertical 
motion analysis. We found that E z of migrators needs to change 
as v m ig(r) ~ exp(-r/2?v), which is the same as the expected 
drop in velocity dispersion. The radial epicyclic frequency, on 
the other hand, decreases as 1 jr. We may, therefore, expect that 
the changes in E z and E p for the same migrating sample are 
different, with migrators contributing more to the radial disc 
heating. 

At this point a deficiency in the simple model we consider 
is found. As revealed by Fig. [5] third row, in our simulations 
outward migrators radially cool the disc inside the bar's OLR 
instead of heating it, i.e., the opposite to the behavior of their 
vertical velocity dispersion. This effect is most prominent in 
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the region between the CR and the OLR for all the models we 
study, while at r > roLR o~ r {r) and cr-(r) behave similarly. 

9.3. The need for a more complex conserved quantity: 
a function of both J z and J p ? 

We found that the conservation of epicyclic vertical action for 
migrators is consistent with our numerical results. In the case 
of the horizontal motion, this simple model is less accurate, as 
we just discussed. This result is not really surprising, as (i) the 
adiabatic approximation is less adequate for horizontal motions 
than vertical ones even in the axisymmetric case 0, and (ii) it 
is also well known that the classical epicycle theory greatly 
underestim ates the radial a ction even in the case of pure planar 
motion (see Dehnenfl 999l for a better approximations of J p ). 



Additionally, the adiabatic invariance of J p is not neces- 
sarily a good approximation for a time evolving system per- 
turbed by multiple asymmetric patterns with stars migrating 
faster than their epicyclic oscillations. We note that the inverse 
correlation between Acr,(r) and A<r r (r) we found in Section|6]is 
suggestive of an exchange between the vertical and radial mo- 
tions. An exchange between J- and J p can be achieved when 
vertical and horizontal resonances coincide in radius. Such lo- 
cations may be found in the region inside the bar's CR and 
between the CR and the OLR, i.e., where we see the largest 



differ ences between Act, and Act,.. For example. ICombes et al 



have shown that the approximate coincidence in radius 



of the bar's vertical and horizontal 2:1 ILRs (inside CR) are re- 
sponsible for the bar's peanut shape. We, therefore, speculate 
that a conserved quantity for migrating stars could be a func- 
tion of both J- and J p , perhaps in a similar fashion to the con- 
servation of the sum J z + J p proposed by ISridhar and Touma 
( 1996bla ) in the context of an adiabatically growing galactic 
disc, which involves a time dependent perturbation as is the 
case of our discs perturbed by multiple patterns^ More work is 
needed to find a possibly similar invariant for the more compli- 
cated case of migrating stars due to disc asymmetries. 



10. Discussion and conclusions 

In this work we have used Tree-SPH N-body (controlled) sim- 
ulations and unconstrained, cosmological re-simulations, to 
study the effect of internal disc evolution on disc thicken- 
ing. All analyzed models have been previously found, or here 
shown, to develop strong multiple patterns, resulting in efficient 
radial disc mixing. 

We devised a procedure for separating migrating from non- 
migrating stars in simulations. We found an excellent corre- 
spondence between the changes in angular momentum of stars 



2 To improve this, iBinnev and McMillanl d201 lh changed the ef- 
fective potential involved in denning the radial epicyclic energy to 
<5>(R,0) + {L z + J-f /R 2 . 

3 Dubbed, "levitation", this mechanism is at work when the 2:2 res- 
onance between the vertical and radial epicyclic frequencies moves 
through phase space, capturing along its path stars with large J p and 
small J z and releasing them into orbits with smaller J r and larger J z 
(and conversely). 



as a function of initial radius and structure in the r^itM ~ r final 
space, thus, verifying that our separating technique is correct 
(see Fig. [4j. Using this simple procedure, we were able, for 
the first time, to study and contrast the kinematics and spatial 
distribution of migrators and non-migrators in galactic discs. 

A number of recent works have assumed that stars migrat- 
ing from the inner disc retain their high energies, thus popu- 
lating a thick disc component. To test this expectation, we first 
compared the disc velocity dispersion profiles in the vertical 
and radial directions, cr.(r) and o>(r), respectively, for the mi- 
grating and non-migrating populations. Surprisingly, we found 
that the migrators' contribution to the vertical velocity disper- 
sion, Acr z , is negative inside the bar's corotation (CR) radius, 
r cr , and less than 5% in the range r a < r < 3hd, where hd is 
the disc's scale-length (see Fig. |6j. This is related to the fact 
that migrating samples change their energy in such a was as to 
match closely the kinematical properties of the non-migrating 
population at the radius of arrival. We thus have shown that, not 
only do migrators not heat significantly, but they cool a large 
portion of the disc. 

Given the migrators' contribution to vertical heating in the 
outer discs and cooling at small radii, disc flaring could be ex- 
pected. We showed that this is indeed the case, by plotting the 
disc scale-height as a function of radius for migrators and non- 
migrators for each simulation. As anticipated, non-migrators 
exhibit an approximately constant scale-height in contrast to 
the flared disc composed of migrators (see Fig. |7). We have to 
note that this flaring is more significant outside the disc's ini- 
tial radial distribution (~ 12 kpc), i.e., where all particles are 
migrators. 

Using a simple galactic disc model, where both the radial, 
J p = E p /k, and vertical, J : = E z /v, actions, are conserved, we 
explained the decrease of velocity dispersion of outward mi- 
grators as disc cooling. In order to conserve the actions of a 
migrating sample of stars, their energies, and thus the velocity 
dispersions, need to change with radius similarly to k and v. 
While stellar samples migrating outwards lose vertical energy, 
(E z ) ~ cr 1 ,, their reduced o% is still higher than what is expected 
from extrapolating cr- ~ 2(r). This is consistent with Fig. [5] 
where we found that the outward migrators slightly heat the 
disc vertically and inward ones cool it. However, instead of a 
larger contribution from migrators to the radial velocity disper- 
sion, as predicted by this simple model, we found that outward 
migrators cool the disc instead of heating it inside about three 
scale-lengths. This means that, in the approximation we con- 
sidered, the radial energy decreases with radius faster than the 
epicyclic frequency, i.e., the cooling is even more drastic than 
predicted by the simple model from Section [9~T| We suggested 
that a better conserved quantity is a function of both the vertical 
and radial actions. 

We estimated that for the gSa and gSO models, bulge contri- 
bution to the thick disc beyond ~ 2.5 disc scale-lengths is only 
about 1%. For a MW-like bulge (the gSb model) no stars reach 
this radius. This is related to the bar's CR lying mostly outside 
such a small bulge (see Fig. [U. It should be kept in mind that 
here the initial distributions of the bulges were prescribed as 
Plummer spheres. As evident form the third row of Fig. a 
strong contribution comes from disc stars inside the bars' CR. 
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Our controlled simulations did not consider gas accretion. 
In a more realistic case, the mass increase due to a continuous 
thin-disc formation would contract the thickened discs and fur- 



ther decrease the scale-heights (IBournaud et a l. , 2009), partic 



ularly in the outer discs, where for the gSb, gSa, and gSO mod- 
els all stars are migrators. To check this, we applied our tech- 
nique to much more realistic, high-resolution, unconstrained 
simulations in the cosmological context, where discs build up 
self-consistently. Interestingly, we found very similar results to 
those obtained with our simpler models, where the migrators' 
contribution to the total increase in random motions is typically 
less then 3%. 

We do not expect insufficient resolution to be affecting our 
conclusions for the reason that low resolution should result in 
more heating, not less, as we find here. Hence, increasing the 
(intermediate) resolution of our controlled simulations would 
give ris e to even less disc thickening Moreover, the recent 
work bv lRoskar et al.l (12012b has shown that increasing the soft- 
ening length from 50 to 500 pc does not result in any signifi- 
cant change in the disc dynamics, therefore, we expect that the 
migration efficiency is not affected. Finally, our results do not 
change significantly when we use the higher resolution cosmo- 
logical re-simulations studied in Sec. [8] In any case, we are not 
interested here so much in the overall heating, but in the differ- 
ence between the increase of velocity dispersion of migrators 
and non-migrators. 

It is remarkable that both our isolated discs and the much 
more sophisticated simulations in the cosmological context 
(which also use a completely different simulation technique) 
exhibit very similar behavior. We believe the inability of migra- 
tion to thicken discs is a fundamental property of secular disc 
evolution, irrespective of the migration mechanism at work. We 
note, however, that when mergers are involved (i.e., when the 
evolution becomes non-secular), the approximate conservation 
of the vertical and radial actions we find here, may not be valid 
anymore. We defer a study of this important topic to a future 
work. 



10.1. Expected signatures of migration in galactic 
discs 

Considering the results of this work, we now outline three pos- 
sible observables for the effects of radial migration. 

(i) Old disc flaring: The effect of migration, driven by in- 
ternal disc evolution, on the disc heating is mostly seen in the 
outer discs. The reason for this is the following. As we found 
in Sec. |H as stars migrate outwards their velocity dispersion 
decreases as to approximately match the energy of stars at the 
final radius (which did not migrate). However, a small amount 
of energy may be retained, with the effect being stronger, the 
farther the stars come from. In other words, while stars do ar- 
rive with slightly larger energy at an outer radius, the effect is 
negligible, unless the change in guiding radius is very large. 
It thus appears that migration would always work as to in- 



4 Note that if the resolution is too high, then the simulations will 
have unrealistic dynamics - real galaxies have sources of noise, such 
as molecular clouds, star clusters, mergers. 



crease/decrease the scale-height in the outer/inner discs, giving 
rise to flaring. We showed in Section|8]that the situation may be 
different when self-consistently growing discs are considered 
(see Fig. [8}. In any case, we do expect that older stellar pop- 
ulations, being exposed to migration mechanisms the longest, 
would exhibit flared vertical profiles, even if the total mass in 
the disc would not flare. Since this effect is related to the con- 
servation of the average vertical action of stars, in the case of 
an active merger history, this signature may have been erased. 

(ii) Flattening of the average vertical a ction radial profile: 
As pointed out by Minchev et al. ( 201 lbl) . an important con- 
sequence of the conservation of actions for migrating stars is 
that radial migration also mixes the radial distributions of the 
actions (J p ) and (/-): the more the disc mixes, the weaker their 
variation with radius. This may offer a way to constrain the 
migration history in the Milky Way. Estimates of the average 
(better conserved) vertical action in observations for different 
populations of stars should reveal different variation with ra- 
dius (flattening) for older groups of stars. 

(iii) Structure in vertical chemical gradients: While we 
here showed that a kinematically thick disc is not expected to 
result from radial migration, signatures of migration may be 
found in the current chemical distribution of stars, since out- 
ward migrators are preferentially deposited at high altitudes, 
the reverse being true for inward ones (see Fig. |5}- Thus, for a 
given chemical evolution model for the Galaxy, predictions can 
be made for vertical chemical gradients as a function of radius 
as observed today. 
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11. Appendix 

7 7.7. Velocity distributions of migrators and 
non-migrators 

In Fig. [3] we showed that as stellar samples migrate, they de- 
crease/increase their velocity dispersions both in the radial and ver- 
tical directions when coming from the inner/outer disc, respectively, 
so as to approximately match the kinematics of non-migrating stars in 
the final radial bin. Fig. [9] offers a different view of this process. We 
consider the gSa model in the time interval 0.5 < t < 1.5 Gyr (same 
as the third column of Fig. [3}- The first row of Fig. |9]shows the radial 
(left) and vertical (right) velocity distributions for stars which stay in 
the radial bin during the entire time. An increase of velocity disper- 
sion is seen with time. The second row shows stars which came from 
r < 6 kpc. Note that their velocity dispersions drop by about 40%. The 
third row plots stars originating from r > 13 kpc. These are found to 
heat up upon arriving at the final radius. In all cases, final distributions 
(orange histogram) roughly match for each velocity component. 

7 1.2. Particle density from selection procedure 

Here we discuss in more detail our procedure for separating migrators 
and non-migrators. Fig.[TT]contrasts the true density (panels 1-3), es- 
timated directly from the initial and final radii of stars in the disc for 
the times shown, and the density obtained in the method described in 
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Section[5](panels 4-6, same as top row of Fig. [4). Note that the overall 
shape of contours remains the same despite the overlap of particles re- 
sulting from our selection method. Considering that the contour levels 
are the same, we see that the true density drops faster with radius. This 
is due to the fact that the overlap of stars involves more stars coming 
from the inner disc, due to the asymmetric drift. It is clear that this 
overlap will not affect velocity measurements since contribution from 
multiple stars would cancel out. In any case, we remove multiples 
prior to estimating the velocity dispersions and disc scale-height. 

11.3. Radial bin overlapping 

To achieve smooth variation with radius when estimating velocity dis- 
persion profiles, etc., for all figures presented in the paper we over- 
lapped bins every 0.3 kpc. Note that a typical radial bin width is about 
2-3 kpc, increasing with radius (see Fig. UOt, In Fig. [13] we examine 
the effect of changing the distance between the radial bins. The top, 
right pair of panels is the same as the top row of Fig. [6] (bin overlap 
of 0.3 kpc), showing the fractional changes in the vertical and radial 
velocity dispersions due to migrators, A<x z (r) and Acr r (r), respectively. 
We contrast four possibilities: bin overlap of 0.1, 0.3, 0.5 and 1 kpc. 
The overall shapes of the functions plotted remain the same, while, 
naturally, structure is lost with an increase in the distance between 
overlapping bins, i.e., for 0.1 to 1 kpc. 

11.4. The effect of different bin sizes 

To construct the figures in this paper we have estimated the bin size as 
2a, where a ~ a, Ik is the approximate epicyclic amplitude at a given 
radius, with cr r and k the radial velocity dispersion and epicyclic fre- 
quency, respectively. In order to contain non-migrating stars, we allow 
for particles to oscillate around their guiding radii by ±a, therefore, we 
select our bin size to be 2a. Since both <x r and k are functions of radius, 
so is the value of a. 

In Fig.[T3]we show how the migrators' contribution to the disc 
velocity dispersion (as shown in Fig. [6} changes with bin size. In ad- 
dition to the choice of 2a used in Fig. [6] (middle two columns here), 
in the first two rows we also consider a bin size of a (columns 1-2) 
and 3a (columns 5-6) for the gSb and gSa models and earlier and later 
times, as indicated. The bottom row shows the effect of a constant bin 
size of 1, 2, and 3 kpc for the gSb model. 

In all cases, the innermost disc appears to be affected the most 
by the change in bin size. Concentrating on the earlier times (when 
most of the velocity dispersion increase takes place), we note that the 
oscillatory inverse functional behavior of Acr z (r) and Acr r (r) is seen 
for all bin sizes. 
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Fig. 9. Distributions of migrators and non-migrators ending up 
in a radial bin in the range 10 < r < 13 kpc, estimated in 
the time period 0.5 < t < 1.5 Gyr for the gSa model. First 
row: The radial (left) and vertical (right) velocity distributions 
for stars which stay in the radial bin during the entire time. 
An increase of velocity dispersion is seen with time. Second 
row: Stars which came from r < 6 kpc. Note that their velocity 
dispersions drop by about 40%. Third row: Stars originting 
from r > 1 3 kpc. These are found to heat up upon arriving at the 
final radius. In all cases, final distributions (orange histogram) 
roughly match for each velocity component. 
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Fig. 13. Changes in Fig.|6]due to the use of different bin sizes in the selection procedure. In addition to our standard choice of 
2a (middle two columns here), in the first two rows we also consider a bin size of a (columns 1-2) and 3a (columns 5-6) for the 
gSb and gSa models and earlier and later times, as indicated. The bottom row shows the effect of a constant bin size of 1, 2, and 
3 kpc for the gSb model. The overall behavior is preserved in all cases. 



